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ABSTRACT 


The illumination condition of an artificial earth satellite 
is investigated in relation to the geographic location of the 
sub-satellite point. A stochastic problem is formulated by 
postulating a "random" satellite, and the solution to this 
problem is indicated analytically. A statistical solution is 
obtained from a digital computer simulation prograrn, together 


with a sensitivity analysis of parametric variations. 
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PREFACE 

The general purpose of this thesis is to investigate certain aspects 
of the illumination of artificial earth satellites by the sun, in relation 
to the track of the sub-satellite point as measured in a geographic 
coordinate system. Imposed upon this study is a stochastic model 
of the earth-satellite system. 

An analytic approach to the basic problem was made, followed 
by a computer simulation program for a statistical solution. The 
basic problem is set forth in the Introduction. 

I wish to express my appreciation to the following persons for 
their contributions to this study: 

Mr. E. J. Leinfelder, Chance-Vought Corp. , Astronautics 
Division, who first posed the problem to me, and who provided the 
stimulus and inspiration leading to this study; 

Professor W. P. Cunningham, U. S. Naval Postgraduate School, 
who indicated the method of formulating the problem in a manner 
which led to the desired results; | 

Professor R. M. Thatcher and Miss Pat Hoang, U. S. Naval 
Postgraduate School, who assisted me in no small measure in the design, 
writing and de-bugging of the computer simulation program; and 

Professor R. R. Read, U. S. Naval Postgraduate School, who 


detected initial errors in the chapter on probability formulation. 








Symbol 


dl 


GLOSSARY 


Definition 
semi-major axis of orbital ellipse 


apparent change in longitude of sub-satellite 
point from initial reference time to time t 


altitude of satellite above surface of earth 


inclination of satellite orbital plane to 
equatorial plane 


inclination of equatorial plane to ecliptic plane 
longitude of designated geographical points 
longitude of sub-satellite point at time t 
longitude of ascending node 

unit vector from center of earth to satellite 
projection of Yr onto equatorial plane 


specified angular distance of sub-satellite point 
from designated geographical points 


instantaneous time after initial reference time 


angular rotational velocity of earth about polar 
axis in radians per unit time 


angular velocity of satellite in radians per unit 
time 


angle between E andr vectors 


angle at satellite subtended by the center of the 
earth and the surface of the earth 


angle between vernal equinox and ascending node 
of satellite 


Vill 








Symbol 


ts | 


1 


M'! 


Definition 


unit vector from center of earth to sun 


angle between ascending node and projection of r 
onto equatorial plane 


angle between ascending node and © vector 
universal gravitational constant 

angle between vernal equinox and E vector 

latitude of designated geographical points 

angle between r andr! 

event that desired passage occurs during first orbit 
mass of the earth 


event that desired illumination occurs at time of 
closest approach to designated geographical point 


radius of the earth 


orbital period of satellite 


Vill 








INTRODUCTION 


It is proposed to investigate some aspects concerning the 
illumination of an artificial earth satellite by the sun in relation to 
the track of the sub-satellite point as measured in a geographic 
coordinate system. But before stating the complete problem, some 
preliminary discussion is necessary. 

Consider a satellite in a circular orbit around the earth. Its 
position may be completely specified at any time (ignoring 
perturbations) if its altitude, orbital inclination to the equatorial 
plane, and the longitude of its ascending node at a particular time is 
known. The longitude of the ascending node may be defined as the 
longitude at which the sub-satellite point crosses the equator entering 
the northern hemisphere. Under these circumstances the track of 
the sub-satellite point on the earth's surface may be computed for any 
set of sequential orbits. 

Due to the rotation of the earth about its polar axis, the longitude 
of the ascending node is displaced westward on each orbit by an 
amount equal to the rotational velocity of the earth times the orbital 
period of the satellite. This displacernent will be denoted w,T, where 
We is the angular rotational velocity of the earth about its polar axis, 
and T is the orbital period of the satellite. 

If it is desired to determine if the track of the sub-~-satellite 


point will pass within a given distance of a given point on the earth's 








surface, then under the stated conditions the problem is completely 
deterministic. A simple but somewhat tedious calculation will 
ascertain the number of orbits required, the exact time, and the 
actual distance from the given point at the time of closest approach. 

Now suppose that, instead of being completely specified, these 
parameters of altitude, inclination, and longitude of ascending node 
were random variables. Conceivably this situation could arise where 
a satellite has been "lost" to surface detection units for an extended 
time period. Similarly, this proposition might be required ina 
feasibility study of a ground visual detection station. In any event, 
under these circumstances the problem of ascertaining if the satel- 
lite track passes within a given distance of a given point is no longer 
deterministic but stochastic. The solution must be given asa 
probability of such passage after a given numberci orbits, and will 
be a function of these three random variables. 

To simplify the problem considerably, let the variables of 
altitude and inclination take on a range of discrete values. Then seek 
a solution for the random variable of longitude of ascending node for 
each pair of values for altitude and inclination. The longitudes of 
the second and al! subsequent ascending nodes are sirmple functions 
of the first nodal longitude, since specifying the altitude specifies 
orbital period, hence specifies the longitude displacement, w,T. 


Then the probability of passage is a function of the random variable 





of longitude of first ascending node, and of the orbital period, T. 

The second part of the problem is to determine if the satellite 
is also illuminated by the sun. It appears obvious that the satellite 
need not be over the sunlit side of the earth to be illuminated. In 
fact, the only requirement for illumination is that the satellite be 
outside the shadow cone of the earth. For satellites at sufficient 
altitudes this may not occur until the satellite is well into the "back 
side" of its orbit. 

With these considerations in mind, the main problem may now 
be stated: 

Given a Satellite at an altitude h, with an orbital inclination i, 
and two geographical points on the earth's surface on the same 
meridian of longitude W, and equidistant from the equator at latitudes 
L North and South: 

Find the joint cumulative probability that the satellite will pass 
within a given angular distance s of either point and that the satellite 
will simultaneously be illuminated by the sun, the result to be stated 
as a cumulative probability distribution function in terms of successive 
orbital periods after initial reference time (arbitrarily set equal to 
zero). 

Some previous studies of satellite illumination have been made, 
made, but no studies have been found which treat the problem froma 


? 


stochastic viewpoint. 








This thesis first develops the necessary mathematical models 
for the orbital mechanics, then sets forth the precise probability 
formulation statements. The statistical solutions are then derived 
from a digital computer simulation program which is discussed in 
the appendices. 

While the design and writing of the computer simulation program 
proved of great interest to the author, and provided an opportunity 
to apply recently acquired knowledge of Monte Carlo techniques, it 
must be emphasized that the computer program is not an entity in 
itself, but merely a useful tool in arriving at a solution to the stated 


problem. 


l 

Leinfelder, E.J., and Thrasher, J.R., Orbit~Earth Shadow 
Relations, AST/EOR-13277, Chance-Vought Corp., Dallas, Texas, 
1960. CONFIDENTIAL 


2Koelle , H.H., Editor, Handbook of Astronautical Engineering, 
McGraw-Hill Book Co. , New York, 1961. 








CHA PVER 


ORBITAL MECHANICS 


Initial As sumptions 


), The earth's solar orbit and the satellite orbit have zero 
eccentricity. 
lar The earth is a perfect sphere with total mass concentrated 


at the center. 

oO The earth's equatorial plane is inclined to the ecliptic 
plane at a constant angle i' (= 23° 26' 59*), 

4. The earth's solar umbra is a circular cylinder (valid for 
values of h less than three times the earth's radius). 

je At initial reference time the longitude of the ascending node 
of the satellite may be considered a random variable, uniformly 
distributed on the interval t 210 | . At this time the satellite will 
be assumed to be located at the ascending node. 

6. All other perturbations of a satellite orbit will be ignored. 


This includes atmospheric drag and solar radiation pressure. 


Discussion of As sumptions 


The effect of these assumptions is to confine the problem to 
one involving only central forces. This also implies that precession 
of the satellite orbital plane will not be included. In order that these 
omissions do not radically affect the solution, the probability dis- 


tribution function sought will only be examined for a pericd of 16 orbits. 





This will be on the order of 24 hours for the range of altitudes con- 
sidered. As will be seen subsequently, two more random variables 
must be introduced; one to account for the space orientation of the 


satellite orbital plane, and another to account for seasonal variations. 


Method of Attack 
A mathematical model will be established representing the motion 
of the satellite and the sun with reference to a geocentric coordinate 


system, followed by a probability formulation. 


Mathematical Model 

Let E represent the vector from the center of the earth to the 
sun, and let r represent the vector from the center of the earth to 
the satellite. The initial assumptions will ensure that these vectors 
are of constant magnitude and possess constant rotational velocities 
about the center of the earth. 

Using Figure I-1, find the direction cosines of E in the rectan- 
gular geocentric coordinate system depicted, where the X-axis is 
along the vernal equinox, the Z axis is directed through the north 
celestial pole, and the Y-axis perpendicular to the other two axes in 
a conventional manner. The origin is, of course, the center of the 
earth. 

During a yearly cycle the angle H takes on all values in the 


interval [0 27 . The values of x, y, and z are easily determined as: 








y 
/\ 
ea roy 


as WEL 
ie a <n 


=n i | eq 
xX 


Fig. [-I- Kepresentation of 
Earth-Sun vector E 








x= Fk, cos H 
y = Ecosi' sinH (1-1) 
z= BE sini' sinH 


Hence the direction cosines are : (cos H, cos i' sin H, sin i' sinH). 

From Figure I[-2 find the direction cosines of r in the same 
coordinate system. The angle D may be thought of as the right 
ascension of the ascending node. The angle G represents the 
angular distance of r from the ascending node. Obviously the angles 
D and G take on all value in lo 21 | ; 

From Figure I-2 the angle F is the projection of G onto the 
equatorial plane. Let r' be the projection of Yr onto the same X-Y 
plane. Then r' = r cos (angle between r andr’). 

Using Figure I-3 the right spherical triangle with sides G, F and 
L '' reveals that L' is the angle between r andr!. Then 

sin L' = sini sinG (1-2) 
and 
cos L' =(1- sin” yee 
Se oeeAc GON (1-3) 
Also from the right spherical triangle in Figure [-3: 
tan F = cositanG ( 1-4) 
hence 


EF = oa (cos i tan G) (1-5) 
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rig.[-2: Kepresentation of Earth 


Satellite vector fF 











Fig. [-5: Simplified view of 
Earth~ Satellite vector ¥ 








Referring again to Figure I-2, 


x= r'cos (D+F) 


ae 


=r cos L' cos (DtF) 


= r (1 - sin i sin’ G) % cos D + ne (eoe i tan a1] ; 
(1-6) 
y = r' sin (D+F) 


=r (]- sae rh Bie Gy? sin [> 42 tan Ycos i tan a1] ; 
(1-7) 


Ze Sin i 


=f Sim a Sin G 


(1-8) 


Using Equations (1-6), (1-7) and (1-8), the direction cosines of r are: 


Ve = 

‘ (1 - aan’ i eine G) * cos D + tan licos i tan <) ; 
Vy = 

(1 - shin” i sin¢ G) * sin |p + tan - (cos i tan G) : 


sini sin G S (1-9) 
Referring to Figure I-4 it is noted that the angle G describes 


the angular position of the satellite (or the sub-satellite point) in the 


orbital plane from the ascending node at time t. 


Hence 
G = w,t (i- 10) 
where 
w, = angular velocity of satellite in radians/uxit time 
t = instantaneous time after reference time. 


As the satellite proceeds in its orbit to timet, the earth also rotates 


about its axis with an angular velocity w, 


Then di is the apparent 


change in the longitude of the sub-satellite point in time t. Clearly, 


itl 
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Fig. 1-4: | ocation of sub- satellite 
point ai) tiunec ae 








F =w,t + dl F (1-11) 


Then the longitude of the sub-satellite point at time t is: 


l 
— 


is + dl 


Totes = wee 
= ee + tan! (cositanG) — w,t (1- 12) 
It is more convenient at this point to develop an analytic model 
which will indicate whether or not the satellite is illuminated by the 
sun. Figure I-5 depicts the plane containing the E andr vectors. 
For simplicity, let both E andr. be unit vectors. 
As the satellite orbits about the earth, the angle A (between 
the E andr vectors) can take on all values in the interval |0 27¢| , 
However, since the representation in Figure I-5 is symmetric about 
the E vector , the problem requires an interest in the magnitude of 
A only in the interval [0 ; rt | 
Clearly, 
E°: r = cosA (1-13) 
If cos A is non-negative, the satellite will be illuminated by the sun. 
If cos A is negative, the satellite may or may not be illuminated. 
Let B = A- 7( /2 as shown in Figure I-5. Let C be the angle at 
the satellite subtended by the center of the earth and the surface of the 


earth. Then 


Ge= sin e€ ’ (1-14) 
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a constant for circular orbits, where ses is the radius oi the earth. 
Then if B + Cis greater than 7(/2, the satellite will be in the 


shadow cylinder of the earth. Hence if the satellite is illuminated, 


then 
oS - 
cos (E° rr) - TH +sin : Bos | = us (1~15) 
Z Re ata 2 
or 
eo aa zs R < 
cos LE * xr) + Sin : e =) oa (1-16) 
Ree tee 


If it is not illuminated, then 


fossa r) + sin?! ae ert (1-17) 


Ie SP 1g) 


Using the results of Equations (1-1) and 1-9), then: 
aie =n | 
jog | ee eee Stn’ Ghee. 
{cos H cos | D a Aen eee i tan G) | 1. 


cos i' sin H sin | D + tan” (cos i tan G) | 


+ sini' sini sin H sinG (1-18) 


The next stage of the problem is to determine if the satellite 
track passes within a given angular distance s of either of two 
geographical points with coordinates (L, W) and (-L,W). With this 
concept the satellite track may be considered to sweep out 2 "swath" 
over the earth's surface, of width 2s. Denote l', as the longitude 


of the ascending node. 


15 








For illustrative purposes it is advantageous to show a schematic 
plot of the satellite track across the surface of the earth, with 
Special emphasis upon its relationship to the designated geographical 
points. For clarity, a Mercator projection will be used. 

For cases where || -s > | , then the situation in 
Figure I-6 will exist. Clearly, if the value of W lies in one of the 
intervals | 1; ; bs forj=1,3, 5, 7, then the cesired passage 
will ecu: Indeed this is a necessary condition. 

For the case where i | ~ Ss = Fa : FH +s, then Figure 
I-7 willapply. Here, the desired passage will occur if W lies in 
one of the intervals [1;. be] , for) = 1, 3. This tooisa 
necessary condition. 

Using these conditions to determine if passage occurs, the time 
t of closest approach to the geographical point may be evaluated. 
Using this value of t to compute cos A in Equation (1-18), then the 
conditions contained in Equation (1-16) are applied to determine if 
the desired illumination occurs. 

It is now necessary to develop expressions for the orbital period 
and angular velocities for the satellite and the earth. It can be shown 

that the period of a satellite is: 


T= 27 (a)?/4 (1-19) 
(amM)jl/2 


16 
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where 


a = semi-major axis of the orbital ellipse 


tI 


Re th tor citeular orbits 

G  =universal gravitational constant 
; 

M =mass of the earth. 


In mks units, 


Eee Ce O° nears 


Re = 
G =6.668 x107!! newton-meter*/kilogram? 
‘ 24. 2; 
M =5.98 X 10°” kilograms 
Then 
Bic 
T = 2 1 (6.371 X 10° meters Raye | 
= \/2 
6.668 x 107!) nt.—m* x 5.98 X 10°* kg. d 
ke 
BW ErA 
= 27 (6.371 x 10°m, +h? 
(3.99 X ioe Sao ee) if 
=37u8 x lon (6.371 X i0° + h) 3/2 SEC 
SN? ee 
= 5 LV (ancien) aoe (1-20) 
where h is given in units of 10° meters. 
Since wot = 2 TT ’ 
we = 27 eal (1~2 1) 
and 
ae 1 Ee 


where T,, is the rotational period of the earth and 


Ji 








[1 


tI 


(2 1/1440) min. 


4.363 X 107° min. 


20 


(1-22) 








CHAPTER II 


PROBABILITY FORMULATION 


The next stage leading to a solution of this problem is the 
formulation of a probability statement. The probability sought is 
the probability that illumination occurs simultaneously with the 


h Orbit, for n= 1) 2. ee 


desired passage by the nt 
Looking at the first orbit, the longitude system may, without any 
loss of generality, be arbitrarily defined such that the longitude of the 


ascending node, l', , is equal to zero, and all other longitudes are 


O 
measured eastward from zero through a total angle of 27T . 

Let M be the event that the desired passage occurs sometime 
during the first orbit, and let N be the event that the desired 
illumination occurs at the particular time t when the satellite makes 
its closest approach to the geographical point. Note that Misa 
time-independent event, while N is strictly dependent upon the value 
of t. The probability sought is P [ oa | 

Studying each event separately, it is clear that 

p[m| =-P[} =w Sere | (2-1) 
where j = 1 and 3, or else j = 1, 3, 5,7, depending upon the values 
ofi, sand L. 

This equation can be written in more meaningful language if 

p | Mj| is defined as the probability that W lies in fi. dar] 


J 


ok 





‘1. <= = _ ai 





Clearly, the events Mj, M3, ..., are mutually exclusive, 
hence 
Pp | M| 7 » Pp | mi, | (eee) 
odd j 
Also, from Equation (1-16), 
< <a 
p[n]=p[a+c S 1 | -P| a = en (2-3) 
Hence the random variable associated with event Mis W, and 
it appears intuitively correct that A is the random variable associated 
with event N. Then the joint probability sought is 


Pp { Mn | - » Pp | MN | (2-4) 


odd j 


To verify that A is a random variable, note that Equation (1-18) 
shows that A is a function of i, G, Dand H. But for a particular 
time t, the values of i and G are constant. Then P gE s Wl -C 
is simply the probability that D and H take on a range of values such 
that A 1 - C, which implies that D and H are random variables 
with associated probability density functions, £p(d) and fyz(h), both 
defined on the interval | 0 oot . Care should be taken to avoid 
confusion with h, the satellite altitude; however, in most Situations 
the meaning will be clear from the context. 

Using the basic assumptions and the physical interpretations of 
D and H, it appears intuitively correct to state that they are inde- 


pendent, uniformly distributed random variables such that 


Ze 






- i) @ sae Me 
sities 
la 
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®@ 


— = —- = —_- & — "aor 





f(a) = | an for0 £a £2a 
(2-5) 
0 elsewhere 
va < 
£,,(h) Ss) | em for0 =h = 2n 
(2-6) 
0 elsewhere 


It quickly follows that A = A(G,H) is a random variable with 
an associated probability density function, f, (2) , defined over the 
interval [0 : TC | : 

While it appears unlikely that the mathematical expression for 
fy (a) could be derived quickly by analytic methods, the formulation 
of the problem will continue as if fq (2) were completely described. 

Before proceeding further with Equation (2-4) it is necessary to 
develop the form of the probability density function for W. Since the 
only quantity of interest is W = l'g , and since one of the basic initial 
assumptions was that l', was uniformly distributed, it follows that 
W is also uniformly distributed. Arbitrarily setting I! equal to 


zero, aS before, merely ensures that W is defined over the interval 


omer | - 


Then 
fywlw) = 1/21 for 0 ewi2nr 
(2-7) 
0 elsewhere 
Now by the rules for conditional probabilities, 
P[ MN | -P|N mM; | P[M,| (2-8) 











and 


p[mn} + > i ive - PLM; | (2-9) 


odd j 
Assume now that, associated with P [fo , there exists a 
conditional probability density function fa \wi4 | w) ,» defined over the 


interval [o : rt | - Then Equation (2-9) may be written as 


1. rw—C 
i= | MN | _ 3 * f En |wl|) f.(w) da dw 


odd j li 0 
fy) TC 
= ] > - | fa |wl2|w) da dw (2-10) 
found 
odd j i 0 


using the results of Equation (2-7). 
The final probability to be calculated is the joint cumulative 


probability that the desired events M and N will occur simultaneously 


th 


orbit, forn=1, 2, ..., 16. To accomplish 


at least once by the n 
this requires some new notation. 
Let P, represent the desired probability of at least one such 


h orbit. Let P(l) be the p=obability that 


joint occurrence by the nt 
the desired event occurs on the first orbit; clearly P({i) is the solution 
to Equation (2-10). Let P(n',nt]) be the joint probability of no such 


occurrence on the first n orbits and one occurrence on the (nt Nee 


orbit. 









» Cc 
on 

ee “im | ces uh 

eel ee nl 

pe ss 6 SS D> | 
= > aa —a- - 

~- ener mey << 

—— a 

=—_> See =e eo 

- 2s Cela 

— ee el 








pag | 











Then 


P, = P()) 
Pe = 2 er(2" 3) 


and, in general, 
Softee ty ode lbaa ate!) ee (2-11) 
This equation describes symbolically the procedure for 
euewlating P. forn=1, 2,..., 16. The problem of caleulatine 
the probabilities for orbits subsequent to the first can be considerably 
alleviated by utilizing a characteristic of the satellite tracks. 
From Equation (1-12) it is easily shown that the longitude of the 
ascending node for the second orbit is 
Vr = Ee = w 0 ; (2~12) 
i,e., I'p has been displaced westward of I', by the amount wT. 
But this is equivalent to retaining l', as the longitude of the ascending 
node for subsequent orbits and displacing W eastward a distance w,T. 
In other words, Equation (2-10) may be used for the nt) orbit, pro- 
viding the longitude of the geographical points is set equal to 
(W), We ta > ewe for n= lnc een (2-13) 
Nonetheless, there remains unsolved the problems of determining 
the limits of integration (with respect to w) and discovering the 


analytic expression for falwilw) in Equation (2-10). These problems 








will be avoided by abandoning the analytic approach in favoy 


of a computer simulation program using Monte Carlo rnethod 


\h 
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CHAPTER Iil 


SIMULATION PROBLEM 


It is expedient first to examine the form of the desired solution 


in order to assess the proper method of arriving at it. 


The cumulative probability distribution function sought is nota 


continuous function, but a monotonic step function, with certain 


probabilities associated with each point corresponding to the end of 


an orbital period. In fact, the desired result will be a distribution 


function truncated at 16 orbital periods. 


Now assume for the moment that the problem can be set up for 


Simulation, where a suitable number, k, of samples can be obtained 


from a pseudo-random number generator. Let 


ie 


= 


a 


me 


Here , 


= number of times that passage occurs on the St orbit 


number of times that passage occurs on the ith orbit 


and not on the i—l previous orbits, for i=2,3,...16. 


= number of times that illumination occurs on the 18t Orbit; 


given that simultaneous passage occurs on the first orbit. 


DD orbit. 


number of times that illumination occurs on the if 
given simultaneous passage on the same orbit, and that 


such joint event did not occur on the i=! previous orbits, 


for i = 223..." 


"passage'' implies that the satellite track passes within the 


rag | 











specified angular distance s of either of the two geographical points. 


With these precise definitions it is possible to assert that: 


E|P,| = —§ [ PCD | = yy 
Ee 
e[P,] =v, + = [Pura} ey, ty, 


and, in general, 


La — » Yj et) 


If Pp) is defined to be the cumulative probability of passage 


(independent of illumination), then, in similar manner, 


a n 
a | Pate) | - 1) a (3-2) 
KLoj=] 


Now define u:; andv:.: such that 


ij ij 
Ui, = l if the desired passage occurs for 
: thy aaa. 
the first time on the i orbit, 
el 
for the j sample random number ? (3-3) 
value | 
0 Otherwise 
Vij = i If the desired illumination an 
passage occurs for the first time om | 


on the ith orbit, for the ,j th cam- » (3-4) 
ple random number value. 


¢) Otherwise 
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where j = 1,2,...,k. Note that v,;, = 1does not require that Uij = i 


J 
Then k 
sist = basi (3-5) 
jet 
k 
and Yi = » Vij (3-6) 
jee 
Hence n k 
E[pP,| 71 . | ‘ (3-7) 
k 
ead] j=l 
and n k 
E [P,(p)] = 2 ) ) Us, (3-8) 
k 
1S) j=l 


Equations (3-7) and G-8) can be considered as the basis for the 
simulation model. 

The design of the computer simulation program is described in 
some detail in Appendix A, with a description of some of the pro- 
cedures used to reduce the required computation time. The language 
used in writing this program was the CDC-1604 version of NELIAC, 
which is peculiarly well-adapted to this type of program. 

Appendix B contains a summary of the output of this program, 
along with a sensitivity analysis for the effect of single parametric 


variations of satellite altitude, orbital inclination, latitude of geo- 


graphical points, and angular distance s. 


29 






7 . 
4 S 





7 s 
; > i * aeeele 
colleen ee oe ee ed 


422 eee aw 
—— ome tt he bones 





BIBLIOGRAPHY 


Space Navigation Handbook, NAVPERS 92988. 


Merrill, G., Editor, Principles of Guided Missile Design, 
D. Van Nostrand Co., Inc., Princeton, New Jersey, 1960. 


Van Allen, J.A., Scientific Uses of Earth Satellites, 
University of Michigan Press, Ann Arbor, Mich. , 1956. 


The NELIAC Cornpiler Language, Department of Operations 


Research, U.S. Naval Postgraduate School, Monterey, Calhf., 
1963. 


30 








APPENDIX A 


DESIGN OF A COMPUTER SIMULATION PROGRAM 


The primary concept of the program is an artificial-time- 
controlled sequence of computations, based upon Equation (1-16) and 
a modification of Figures I-6 and I-7. The satellite is figuratively 
put in orbit at its ascending node, corresponding to time t = 0, and 
where the longitude 1', is also identically zero. A suitable time 
increment is used to activate an artificial time "clock" which 
determines the instantaneous latitude and longitude of the sub- 
satellite point. 

Certain short~cuts are immediately available to reduce the 
otherwise inordinate number of calculations required. 

Figure A-lis similar to Figure I-6, except that the end points of 
the longitude intervals correspond to the instantaneous longitudes of 
the satellite track when the instantaneous latitudes equal either L-s 
or L+s. Clearly the desired passage cannot occur unless the 
instantaneous track latitude L'; lies in the interval [L-s : Lts| : 
which implies that if the instantaneous longitude l', lies in the interval 

E , 141 for j= 1, 3, 5, 7, then passage will occur if W also lies 
within this same interval. 

In similar manner, the interval EF 143 ign 1 = 1,5, is defined 


rn“ . : : 
for situations where L-s =i = Lts. See Figure A-2. Again, 
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clearly, passage will occur if l'; and W lie within this interval. The 
numbering of the subscript j is designed to facilitate the computer 
programming. 


From Equations (1-2) and (1-10): 


sin L' = sini sin wt 
hence for L' = Legs, 
t = 1 sin | sin (L —s) 7 (A -1) 
ie sini 


and for L't = Lts, 


1 sin sin (L + s) =to (A -2) 
Ww sini 


where t; corresponds to the time at which the satellite track attains 


J 


a value of longitude equal to li}. The remaining values of t; for j 


2 can be determined quickly by 


a Ji = to 
Z 

t, = Yes ty 
Z 

ioe kiueeteents 
2 

[aoe St it 
Zz 

t2 = T + tft, 
Z 
a 


which can be easily verified by the symmetry of the track. T is of 
course the satellite period. 


Using Equation (1-12) it follows that 


A-4 


| 
L = tan (cos i tan W tj) - vel, 


and accordingly 
eee tan7\cos itan wet) -w,t] (A -3) 
I= tan "cos i tan Wet2) a w to (A -4) 
where l,and 1, are defined as in Figures A-land A-2. In these 
figures if W lies within any of the intervals depicted then passage 
will occur. The remaining values of i can be computed in exactly 
the same manner as Equations (A-3) and (A-4). 

Then the first test for the desired passage is determining if the 
sample value for the random variable W lies within one of the afore- 
mentioned intervals. If W satisfies this test, then the test for illumi- 
nation may be made. However, it is possible for W to be outside 
such intervals and still provide for the desired passage. If W lies 
outside of, but sufficiently close to, either end point of the depicted 
intervals, then passage may occur as in Figure A-3. 

At this point it is appropriate to introduce the approximation 
s ' of the true instantaneous angular distance between the sub-satellite 
point and the designated geographical point. For values of s' 
sufficiently small, the curved surface of the earth can be approximated 
by a plane. Then, from Figure A-3, 


“| 


1/2 
L Li) 7 + |w -1,|? cos? | (A -5) 








———— ~_ —4e 








Fig.A-3: Mercator projection ot track: 
special case 








where the cosine term is used to compensate for the "shrinkage"! 
in distance between adjacent meridians of longitude with increasing 
latitude. 

It should be obvious that passage cannot occur if |(W -l';) x 

mee 
cos L;| =s, hence this quantity may be used as an upper bound 
in an iterative procedure for computing successive values of s' for 
successive values of time. Due to the monotonic nature of the 
satellite track, a minimum value of s' will exist for some value oft, 
; : < 

say t*, This value of s' is compared to s; if s' =s, then the desired 


passage will occur. Otherwise, passage cannot occur. 


h th 


If these tests are applied to the im orbit and j~° sample value of 
the random variable W, the value of Uij in Equation (3-3) is determined. 
If uij = 1, then the illumination condition may be determined by 
application of Equation (1-16). Applying this test to the ith orbit and 
Me sample values of the random variables W, D and H determines 
the value of vj; in Equation (3-4). 

However, it may happen that the desired passage will occur with 
no illumination on an earlier orbit, but with the desired passage and 
illumination on a later orbit. In order to compute P,,(p), it is 
necessary to record the earlier event where ujj = 1, but continuing 
the process in order to compute P,. 


Since this simulation program was only incidental to the solution 


of the main problem, the computer flow diagrams and NELIAC 





language flow charts have not been included. However, the IBM 
binary card decks and program print-out are available from 
Professor Richard M. Thatcher, of the Operations Analysis 


Department. 








APPENDIX B 


SUMMARY OF PROGRAM OUTPUT 


The parametric inputs to the simulation program were satellite 
altitude h, orbital inclination i, latitude of designated geographical 
points L, and angular distance from these points s. The values 
of h and s are given in units of nautical miles, but were converted 
internally to units of 10° meters and radians respectively; the 
values of i and L are given in degrees, but were likewise converted 
internally to radians. The values of these parameters and the in- 


cremental values used in the program are given below: 


h: 100(300)1900 nautical miles 
i: 20°(8°)68° 

L: 25°(10°)65° 

S: 100(100)500 nautical miles 


The value of k (sample size) was hopefully selected to be 200, 
in order that the approximations contained in Chapter III might be 
statistically valid, and yet require a reasonably short amonnt of 
computer time. Even so, about three and one-half hours actual 
operating time was consumed on the computer for a single run 
through all values of the above parameters. 

One immediate result from the output of this program was the 


discovery of a "semi-optimal" distance s = 300 nautical miles. 
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In many instances, for various values of h, i and L, the cumulative 
probabilities revealed little or no increase with values of s greater 
than 300. Ina significant number of cases the rate of increase in 
probability declined sharply after s = 300; it is to be noted, however, 
that this was not always the case. Frequently the rate of increase 
in probability with increasing s was approximately linear. 
Nevertheless, it was decided to limit the tabulated values of the 
probabilities to those computed for s = 300. To have done otherwise 
would have increased the size of this appendix by a factor of five. 

Another unexpected discovery concerned certain substantial 
decreases in cumulative probabilities for h = 1900 nautical miles. 
The reduction did not occur consistently for all values of the other 
parameters, however. It is suspected that this reduction is most 
likely a statistical anomaly; unfortunately insufficient time was 
available for further repeated runs uSing this value of h. This so- 
called anomaly is suggested as a topic for further study, as is an 
analysis of variance for repeated runs. The computer program is 
easily adapted to the study of the effects of variation of any one 
parameter. 

Figures B-lthrough B-10 give the joint probability of passage 
and illumination vs. number of orbital periods, while Figures 
B-1l]l through B-15 give the independent cumulative probability of 


passage vs. number of orbital periods. 
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Figures B-|l through B-4 and B-ll through B-15 represent 
probabilities averaged for all the altitudes examined and, as such, 
show the sensitivity due to variations in inclination. The rationale 
for this averaging process lies in the fact that the variations in 
probabilities due to increasing altitude were sufficiently small to 
warrant doing so. 

Figures B-5 through B-1l0, however, are given for specific 
latitudes and inclinations, and show the sensitivities due to 
variations in altitude. In these cases it will be noted that the values 
of Landi are approximately the same. In these cases the probabi- 
lities are consistently higher than in cases where i is considerably 
greater than L. Of course, wheniis less than L by 3° or more, 
no passage occurs. 

These graphs have been smoothed to some extent over the raw 
statistical output, but the general trend has been retained in all 
cases. No attempt should be made to interpolate between discrete 


orbital periods. 
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